<!DOCTYPE html>
<html lang='zh-CN'>

<head>
  <meta name="generator" content="Hexo 5.4.0">
  <meta charset="utf-8">
  

  <meta http-equiv='x-dns-prefetch-control' content='on' />
  <link rel='dns-prefetch' href='https://cdn.jsdelivr.net'>
  <link rel="preconnect" href="https://cdn.jsdelivr.net" crossorigin>
  <link href="https://cdn.bootcdn.net/ajax/libs/KaTeX/0.13.0/katex.min.css" rel="stylesheet">
	<script src="https://cdn.bootcdn.net/ajax/libs/KaTeX/0.13.0/katex.min.js"></script>
	<script>
    
		var katex_config = {
	  delimiters: 
	  [
		  {left: "$$", right: "$$", display: true},
  		{left: "$", right: "$", display: false}
	  ]
    };
    


	</script>
	<script defer src="https://cdn.bootcdn.net/ajax/libs/KaTeX/0.13.0/contrib/auto-render.min.js" onload="renderMathInElement(document.body,katex_config)"></script>

  <meta name="renderer" content="webkit">
  <meta name="force-rendering" content="webkit">
  <meta http-equiv="X-UA-Compatible" content="IE=Edge,chrome=1">
  <meta name="HandheldFriendly" content="True" >
  <meta name="apple-mobile-web-app-capable" content="yes">
  <meta name="viewport" content="width=device-width, initial-scale=1, maximum-scale=1">

  <title>使用 Python 实现素数筛时发现的有趣问题 - Debug Myself</title>

  

  
    <meta name="description" content="最近学习素数筛算法，算法参考自两篇 CSDN 的文章12，原文是用 C&#x2F;C++ 来写的，我用 Python 进行了改写。 12345import timefrom icecream import icimport numbafrom numba import jitimport numpy as np 这里我使用了 icecream 模块进行调试。 先简单介绍一下素数筛算法。对于一个给定的自然数序">
<meta property="og:type" content="article">
<meta property="og:title" content="使用 Python 实现素数筛时发现的有趣问题">
<meta property="og:url" content="https://stassenger.gitee.io/2021/07/08/a%20interesting%20problem%20when%20I%E2%80%98m%20learning%20Euler%20and%20Eratosthenes%20method/index.html">
<meta property="og:site_name" content="Debug Myself">
<meta property="og:description" content="最近学习素数筛算法，算法参考自两篇 CSDN 的文章12，原文是用 C&#x2F;C++ 来写的，我用 Python 进行了改写。 12345import timefrom icecream import icimport numbafrom numba import jitimport numpy as np 这里我使用了 icecream 模块进行调试。 先简单介绍一下素数筛算法。对于一个给定的自然数序">
<meta property="og:locale" content="zh_CN">
<meta property="og:image" content="https://stassenger.gitee.io/2021/07/08/a%20interesting%20problem%20when%20I%E2%80%98m%20learning%20Euler%20and%20Eratosthenes%20method/nojit.svg">
<meta property="og:image" content="https://stassenger.gitee.io/2021/07/08/a%20interesting%20problem%20when%20I%E2%80%98m%20learning%20Euler%20and%20Eratosthenes%20method/jit.svg">
<meta property="article:published_time" content="2021-07-08T05:14:00.000Z">
<meta property="article:modified_time" content="2021-08-14T02:28:03.479Z">
<meta property="article:author" content="Asante Stassenger">
<meta property="article:tag" content="Python">
<meta property="article:tag" content="Algorithm">
<meta name="twitter:card" content="summary">
<meta name="twitter:image" content="https://stassenger.gitee.io/2021/07/08/a%20interesting%20problem%20when%20I%E2%80%98m%20learning%20Euler%20and%20Eratosthenes%20method/nojit.svg">
  
  

  <!-- feed -->
  

  
    
<link rel="stylesheet" href="/css/main.css">

  

  
<!-- hexo injector head_end start -->
<link rel="stylesheet" href="https://cdn.jsdelivr.net/npm/katex@0.12.0/dist/katex.min.css">

<link rel="stylesheet" href="https://cdn.jsdelivr.net/npm/hexo-math@4.0.0/dist/style.css">
<!-- hexo injector head_end end --></head>

<body>
  


  <div class='l_body' id='start'>
    <aside class='l_left'>
    

<header class='header'>
  <div class='logo-wrap'>
  
  
    <a class='title' href='/'>
      Debug Myself
    </a>
  
</div>

  <nav class="menu dis-select"><a class="nav-item active" href="/">文章</a><a class="nav-item" href="/2021/08/09/about/">关于</a></nav>
</header>

<div class='widgets'>
  
    
      
      

<div class="widget-wrap" id="recent"><div class="widget-header h4 dis-select"><span class="name">最近更新</span></div><div class="widget-body fs14"><div class="post-title reveal"><div class="cap"><time>2021-08-14</time></div><a href="/2021/07/08/a%20interesting%20problem%20when%20I%E2%80%98m%20learning%20Euler%20and%20Eratosthenes%20method/">使用 Python 实现素数筛时发现的有趣问题</a></div><div class="post-title reveal"><div class="cap"><time>2021-08-14</time></div><a href="/2021/08/08/how-to-automatically-report-nju-health-information-on-the-phone/">如何在手机上设置南大 APP 自动打卡（一）</a></div><div class="post-title reveal"><div class="cap"><time>2021-08-09</time></div><a href="/2021/08/09/about/">关于</a></div><div class="post-title reveal"><div class="cap"><time>2021-08-09</time></div><a href="/2021/06/12/a%20problem%20when%20adding%20NJU%20email%20to%20Window%2010%20mail%20app/">关于在 Windows10 自带的邮件应用中添加南大邮箱遇到的问题</a></div><div class="post-title reveal"><div class="cap"><time>2021-07-12</time></div><a href="/2021/07/12/sort-algorithm-reminder-3/">排序算法备忘录（三）</a></div></div></div>


    
  
</div>


    </aside>
    <div class='l_main'>
      

      


  <div class='bread-nav fs12'>
  
    
    <div id='breadcrumb'>
      <a class='cap breadcrumb' href='/'>主页</a>
      <span class="sep"></span>
      <a class='cap breadcrumb' href='/'>文章</a>
      
    </div>
    <div id='post-meta'>
      发布于&nbsp;<time datetime='2021-07-08T05:14:00.000Z'>2021-07-08</time>
    </div>
  
  </div>


<article class='content md post'>
<h1 class="article-title"><span>使用 Python 实现素数筛时发现的有趣问题</span></h1>
<p>最近学习素数筛算法，算法参考自两篇 CSDN 的文章<sup id="fnref:1"><a href="#fn:1" rel="footnote">1</a></sup><sup id="fnref:2"><a href="#fn:2" rel="footnote">2</a></sup>，原文是用 C/C++ 来写的，我用 Python 进行了改写。</p>
<figure class="highlight python"><table><tr><td class="gutter"><pre><span class="line">1</span><br><span class="line">2</span><br><span class="line">3</span><br><span class="line">4</span><br><span class="line">5</span><br></pre></td><td class="code"><pre><span class="line"><span class="keyword">import</span> time</span><br><span class="line"><span class="keyword">from</span> icecream <span class="keyword">import</span> ic</span><br><span class="line"><span class="keyword">import</span> numba</span><br><span class="line"><span class="keyword">from</span> numba <span class="keyword">import</span> jit</span><br><span class="line"><span class="keyword">import</span> numpy <span class="keyword">as</span> np</span><br></pre></td></tr></table></figure>
<p>这里我使用了 <code>icecream</code> 模块进行调试。</p>
<p>先简单介绍一下素数筛算法。对于一个给定的自然数序列，要找出其中的所有素数，可以通过筛除合数的方法来寻找。筛法中思路最简单的是埃氏（Eratosthenes）筛，基本思路是：若要寻找 n 以内的所有素数，只需找不大于 <span class="katex"><span class="katex-mathml"><math><semantics><mrow><msup><mi>n</mi><mrow><mn>1</mn><mi mathvariant="normal">/</mi><mn>2</mn></mrow></msup></mrow><annotation encoding="application/x-tex">n^{1/2}</annotation></semantics></math></span><span class="katex-html" aria-hidden="true"><span class="base"><span class="strut" style="height:0.8879999999999999em;vertical-align:0em;"></span><span class="mord"><span class="mord mathdefault">n</span><span class="msupsub"><span class="vlist-t"><span class="vlist-r"><span class="vlist" style="height:0.8879999999999999em;"><span style="top:-3.063em;margin-right:0.05em;"><span class="pstrut" style="height:2.7em;"></span><span class="sizing reset-size6 size3 mtight"><span class="mord mtight"><span class="mord mtight">1</span><span class="mord mtight">/</span><span class="mord mtight">2</span></span></span></span></span></span></span></span></span></span></span></span> 的素数，剔除它们的 2 倍或 2 倍以上的数，从小到大，依次剔除，即可不用单独判断素数，得到素数表。但是这种算法存在一个缺陷，即很多数的质因数不止一个，所以它们被重复判断了。于是有了改进后的线性筛，即欧拉（Euler）筛。欧拉筛的原理是对于每次筛除，只用它的最小质因数筛去，即把这个数表示为一个较大数的质数倍数，这样每个数只用筛去一次，消除了重复。</p>
<p>以下是我从 C 改写的 Python 代码：</p>
<figure class="highlight python"><table><tr><td class="gutter"><pre><span class="line">1</span><br><span class="line">2</span><br><span class="line">3</span><br><span class="line">4</span><br><span class="line">5</span><br><span class="line">6</span><br><span class="line">7</span><br><span class="line">8</span><br><span class="line">9</span><br><span class="line">10</span><br><span class="line">11</span><br><span class="line">12</span><br><span class="line">13</span><br><span class="line">14</span><br><span class="line">15</span><br><span class="line">16</span><br><span class="line">17</span><br><span class="line">18</span><br><span class="line">19</span><br><span class="line">20</span><br><span class="line">21</span><br><span class="line">22</span><br><span class="line">23</span><br><span class="line">24</span><br><span class="line">25</span><br><span class="line">26</span><br><span class="line">27</span><br><span class="line">28</span><br><span class="line">29</span><br><span class="line">30</span><br><span class="line">31</span><br><span class="line">32</span><br><span class="line">33</span><br><span class="line">34</span><br></pre></td><td class="code"><pre><span class="line"><span class="function"><span class="keyword">def</span> <span class="title">Eratosthenes</span>(<span class="params">n</span>):</span></span><br><span class="line">    primes = [<span class="number">1</span> <span class="keyword">for</span> _ <span class="keyword">in</span> <span class="built_in">range</span>(n+<span class="number">1</span>)]</span><br><span class="line">    primes[<span class="number">0</span>] = primes[<span class="number">1</span>] = <span class="number">0</span></span><br><span class="line">    count = <span class="number">0</span></span><br><span class="line">    <span class="keyword">for</span> i <span class="keyword">in</span> <span class="built_in">range</span>(<span class="number">2</span>, n + <span class="number">1</span>):</span><br><span class="line">        <span class="keyword">if</span> primes[i]:</span><br><span class="line">            <span class="keyword">for</span> j <span class="keyword">in</span> <span class="built_in">range</span>(i * <span class="number">2</span>, n + <span class="number">1</span>, i):</span><br><span class="line">                primes[j] = <span class="number">0</span></span><br><span class="line">                count += <span class="number">1</span></span><br><span class="line">    c = <span class="number">0</span></span><br><span class="line">    <span class="keyword">for</span> i <span class="keyword">in</span> primes:</span><br><span class="line">        <span class="keyword">if</span> i:</span><br><span class="line">            c += <span class="number">1</span></span><br><span class="line">    <span class="keyword">return</span> count,c</span><br><span class="line">  </span><br><span class="line"><span class="function"><span class="keyword">def</span> <span class="title">Euler</span>(<span class="params">n</span>):</span></span><br><span class="line">    primes = [<span class="number">1</span> <span class="keyword">for</span> _ <span class="keyword">in</span> <span class="built_in">range</span>(n+<span class="number">1</span>)]</span><br><span class="line">    primes[<span class="number">0</span>] = primes[<span class="number">1</span>] = <span class="number">0</span></span><br><span class="line">    res = []</span><br><span class="line">    count = <span class="number">0</span></span><br><span class="line">    <span class="keyword">for</span> i <span class="keyword">in</span> <span class="built_in">range</span>(<span class="number">2</span>, n + <span class="number">1</span>):</span><br><span class="line">        <span class="keyword">if</span> primes[i]:</span><br><span class="line">            res.append(i)     </span><br><span class="line">        j = <span class="number">0</span></span><br><span class="line">        index = res[j] * i</span><br><span class="line">        <span class="keyword">while</span> j &lt; <span class="built_in">len</span>(res) <span class="keyword">and</span> index &lt;= n:</span><br><span class="line">            count += <span class="number">1</span></span><br><span class="line">            primes[index] = <span class="number">0</span></span><br><span class="line">            <span class="keyword">if</span> i % res[j] == <span class="number">0</span>:</span><br><span class="line">                <span class="keyword">break</span></span><br><span class="line">            <span class="keyword">else</span>:</span><br><span class="line">                j += <span class="number">1</span></span><br><span class="line">                index = res[j] * i</span><br><span class="line">    <span class="keyword">return</span> count,<span class="built_in">len</span>(res) </span><br></pre></td></tr></table></figure>
<p>其中第一个返回值是筛除次数，第二个返回值是素数个数。</p>
<p>由于 Python 中数组的处理可以使用 Numpy，所以我又对代码进行了改进，得到了一系列函数：</p>
<p>首先，对于 Eratosthenes 方法，从 <span class="katex"><span class="katex-mathml"><math><semantics><mrow><mn>2</mn><mi>i</mi></mrow><annotation encoding="application/x-tex">2i</annotation></semantics></math></span><span class="katex-html" aria-hidden="true"><span class="base"><span class="strut" style="height:0.65952em;vertical-align:0em;"></span><span class="mord">2</span><span class="mord mathdefault">i</span></span></span></span> 开始显然是低效的，改为从 <span class="katex"><span class="katex-mathml"><math><semantics><mrow><msup><mi>i</mi><mn>2</mn></msup></mrow><annotation encoding="application/x-tex">i^2</annotation></semantics></math></span><span class="katex-html" aria-hidden="true"><span class="base"><span class="strut" style="height:0.8141079999999999em;vertical-align:0em;"></span><span class="mord"><span class="mord mathdefault">i</span><span class="msupsub"><span class="vlist-t"><span class="vlist-r"><span class="vlist" style="height:0.8141079999999999em;"><span style="top:-3.063em;margin-right:0.05em;"><span class="pstrut" style="height:2.7em;"></span><span class="sizing reset-size6 size3 mtight"><span class="mord mtight">2</span></span></span></span></span></span></span></span></span></span></span> 开始：</p>
<figure class="highlight python"><table><tr><td class="gutter"><pre><span class="line">1</span><br><span class="line">2</span><br><span class="line">3</span><br><span class="line">4</span><br><span class="line">5</span><br><span class="line">6</span><br><span class="line">7</span><br><span class="line">8</span><br><span class="line">9</span><br><span class="line">10</span><br><span class="line">11</span><br><span class="line">12</span><br><span class="line">13</span><br><span class="line">14</span><br></pre></td><td class="code"><pre><span class="line"><span class="function"><span class="keyword">def</span> <span class="title">Eratosthenes_improved</span>(<span class="params">n</span>):</span></span><br><span class="line">    primes = [<span class="number">1</span> <span class="keyword">for</span> _ <span class="keyword">in</span> <span class="built_in">range</span>(n+<span class="number">1</span>)]</span><br><span class="line">    primes[<span class="number">0</span>] = primes[<span class="number">1</span>] = <span class="number">0</span></span><br><span class="line">    count = <span class="number">0</span></span><br><span class="line">    <span class="keyword">for</span> i <span class="keyword">in</span> <span class="built_in">range</span>(<span class="number">2</span>, n + <span class="number">1</span>):</span><br><span class="line">        <span class="keyword">if</span> primes[i]:</span><br><span class="line">            <span class="keyword">for</span> j <span class="keyword">in</span> <span class="built_in">range</span>(i * i, n + <span class="number">1</span>, i):</span><br><span class="line">                primes[j] = <span class="number">0</span></span><br><span class="line">                count += <span class="number">1</span></span><br><span class="line">    c = <span class="number">0</span></span><br><span class="line">    <span class="keyword">for</span> i <span class="keyword">in</span> primes:</span><br><span class="line">        <span class="keyword">if</span> i:</span><br><span class="line">            c += <span class="number">1</span></span><br><span class="line">    <span class="keyword">return</span> count,c</span><br></pre></td></tr></table></figure>
<p>考虑使用 Numpy。事实上，之后的实验发现，不恰当的使用 Numpy 并不会提升速度。</p>
<figure class="highlight python"><table><tr><td class="gutter"><pre><span class="line">1</span><br><span class="line">2</span><br><span class="line">3</span><br><span class="line">4</span><br><span class="line">5</span><br><span class="line">6</span><br><span class="line">7</span><br><span class="line">8</span><br><span class="line">9</span><br><span class="line">10</span><br><span class="line">11</span><br><span class="line">12</span><br><span class="line">13</span><br><span class="line">14</span><br></pre></td><td class="code"><pre><span class="line"><span class="function"><span class="keyword">def</span> <span class="title">Eratosthenes_numpy</span>(<span class="params">n</span>):</span></span><br><span class="line">    primes = np.ones(n + <span class="number">1</span>)</span><br><span class="line">    primes[<span class="number">0</span>] = primes[<span class="number">1</span>] = <span class="number">0</span></span><br><span class="line">    count = <span class="number">0</span></span><br><span class="line">    <span class="keyword">for</span> i <span class="keyword">in</span> <span class="built_in">range</span>(<span class="number">2</span>, <span class="built_in">int</span>(n ** <span class="number">0.5</span>) + <span class="number">1</span>):</span><br><span class="line">        <span class="keyword">if</span> primes[i]:</span><br><span class="line">            <span class="keyword">for</span> j <span class="keyword">in</span> <span class="built_in">range</span>(i * i, n + <span class="number">1</span>, i):</span><br><span class="line">                primes[j] = <span class="number">0</span></span><br><span class="line">                count += <span class="number">1</span></span><br><span class="line">    c = <span class="number">0</span></span><br><span class="line">    <span class="keyword">for</span> i <span class="keyword">in</span> primes:</span><br><span class="line">        <span class="keyword">if</span> i:</span><br><span class="line">            c += <span class="number">1</span></span><br><span class="line">    <span class="keyword">return</span> count,c</span><br></pre></td></tr></table></figure>
<p>对于 Euler 方法，append 的使用会影响速度。事实上，在埃氏筛的算法中，缺少了生成素数列表的操作，额外的实验表明 ，这些操作所占的时间并不会影响整体的结果。不过，这里还是对 <code>append</code> 进行了改进，另外也消除了 <code>len</code> 对于速度的影响：</p>
<figure class="highlight python"><table><tr><td class="gutter"><pre><span class="line">1</span><br><span class="line">2</span><br><span class="line">3</span><br><span class="line">4</span><br><span class="line">5</span><br><span class="line">6</span><br><span class="line">7</span><br><span class="line">8</span><br><span class="line">9</span><br><span class="line">10</span><br><span class="line">11</span><br><span class="line">12</span><br><span class="line">13</span><br><span class="line">14</span><br><span class="line">15</span><br><span class="line">16</span><br><span class="line">17</span><br><span class="line">18</span><br><span class="line">19</span><br><span class="line">20</span><br><span class="line">21</span><br></pre></td><td class="code"><pre><span class="line"><span class="function"><span class="keyword">def</span> <span class="title">Euler_improved</span>(<span class="params">n</span>):</span></span><br><span class="line">    primes = [<span class="number">1</span> <span class="keyword">for</span> _ <span class="keyword">in</span> <span class="built_in">range</span>(n + <span class="number">1</span>)]</span><br><span class="line">    primes[<span class="number">0</span>] = primes[<span class="number">1</span>] = <span class="number">0</span></span><br><span class="line">    res = [<span class="number">1</span> <span class="keyword">for</span> _ <span class="keyword">in</span> <span class="built_in">range</span>(<span class="built_in">int</span>(n / <span class="number">2</span>))]</span><br><span class="line">    count = <span class="number">0</span></span><br><span class="line">    c = <span class="number">0</span></span><br><span class="line">    <span class="keyword">for</span> i <span class="keyword">in</span> <span class="built_in">range</span>(<span class="number">2</span>, n + <span class="number">1</span>):</span><br><span class="line">        <span class="keyword">if</span> primes[i]:</span><br><span class="line">            res[c] = i</span><br><span class="line">            c += <span class="number">1</span>  </span><br><span class="line">        j = <span class="number">0</span></span><br><span class="line">        index = res[j] * i</span><br><span class="line">        <span class="keyword">while</span> j &lt; c <span class="keyword">and</span> index &lt;= n:</span><br><span class="line">            count += <span class="number">1</span></span><br><span class="line">            primes[index] = <span class="number">0</span></span><br><span class="line">            <span class="keyword">if</span> i % res[j] == <span class="number">0</span>:</span><br><span class="line">                <span class="keyword">break</span></span><br><span class="line">            <span class="keyword">else</span>:</span><br><span class="line">                j += <span class="number">1</span></span><br><span class="line">                index = res[j] * i</span><br><span class="line">    <span class="keyword">return</span> count</span><br></pre></td></tr></table></figure>
<p>考虑使用 Numpy 的情况，首先是保留<code>append</code>，之后的实验表明，Numba 的加速会产生相反的结果：</p>
<figure class="highlight python"><table><tr><td class="gutter"><pre><span class="line">1</span><br><span class="line">2</span><br><span class="line">3</span><br><span class="line">4</span><br><span class="line">5</span><br><span class="line">6</span><br><span class="line">7</span><br><span class="line">8</span><br><span class="line">9</span><br><span class="line">10</span><br><span class="line">11</span><br><span class="line">12</span><br><span class="line">13</span><br><span class="line">14</span><br><span class="line">15</span><br><span class="line">16</span><br><span class="line">17</span><br><span class="line">18</span><br><span class="line">19</span><br><span class="line">20</span><br></pre></td><td class="code"><pre><span class="line"><span class="function"><span class="keyword">def</span> <span class="title">Euler_numpy</span>(<span class="params">n</span>):</span></span><br><span class="line">    primes = np.ones(n + <span class="number">1</span>)</span><br><span class="line">    primes[<span class="number">0</span>] = primes[<span class="number">1</span>] = <span class="number">0</span></span><br><span class="line">    res = []</span><br><span class="line">    count = <span class="number">0</span></span><br><span class="line"></span><br><span class="line">    <span class="keyword">for</span> i <span class="keyword">in</span> <span class="built_in">range</span>(<span class="number">2</span>, n + <span class="number">1</span>):</span><br><span class="line">        <span class="keyword">if</span> primes[i]:</span><br><span class="line">            res.append(i)</span><br><span class="line">        j = <span class="number">0</span></span><br><span class="line">        index = res[j] * i</span><br><span class="line">        <span class="keyword">while</span> j &lt; <span class="built_in">len</span>(res) <span class="keyword">and</span> index &lt;= n:</span><br><span class="line">            count += <span class="number">1</span></span><br><span class="line">            primes[index] = <span class="number">0</span></span><br><span class="line">            <span class="keyword">if</span> np.mod(i, res[j]) == <span class="number">0</span>:</span><br><span class="line">                <span class="keyword">break</span></span><br><span class="line">            <span class="keyword">else</span>:</span><br><span class="line">                j += <span class="number">1</span></span><br><span class="line">                index = res[j] * i</span><br><span class="line">    <span class="keyword">return</span> count</span><br></pre></td></tr></table></figure>
<p>接下来是改进的情况：</p>
<figure class="highlight python"><table><tr><td class="gutter"><pre><span class="line">1</span><br><span class="line">2</span><br><span class="line">3</span><br><span class="line">4</span><br><span class="line">5</span><br><span class="line">6</span><br><span class="line">7</span><br><span class="line">8</span><br><span class="line">9</span><br><span class="line">10</span><br><span class="line">11</span><br><span class="line">12</span><br><span class="line">13</span><br><span class="line">14</span><br><span class="line">15</span><br><span class="line">16</span><br><span class="line">17</span><br><span class="line">18</span><br><span class="line">19</span><br><span class="line">20</span><br><span class="line">21</span><br><span class="line">22</span><br></pre></td><td class="code"><pre><span class="line"><span class="function"><span class="keyword">def</span> <span class="title">Euler_improved_numpy</span>(<span class="params">n</span>):</span></span><br><span class="line">    primes = np.ones(n + <span class="number">1</span>)</span><br><span class="line">    primes[<span class="number">0</span>] = primes[<span class="number">1</span>] = <span class="number">0</span></span><br><span class="line">    res = np.zeros(<span class="built_in">int</span>(n / <span class="number">2</span>))</span><br><span class="line">    count = <span class="number">0</span></span><br><span class="line">    c = <span class="number">0</span></span><br><span class="line">    <span class="keyword">for</span> i <span class="keyword">in</span> <span class="built_in">range</span>(<span class="number">2</span>, n + <span class="number">1</span>):</span><br><span class="line">        <span class="keyword">if</span> primes[i]:</span><br><span class="line">            res[c] = i</span><br><span class="line">            c += <span class="number">1</span></span><br><span class="line">        </span><br><span class="line">        j = <span class="number">0</span></span><br><span class="line">        index = <span class="built_in">int</span>(res[j] * i)</span><br><span class="line">        <span class="keyword">while</span> j &lt; c <span class="keyword">and</span> index &lt;= n:</span><br><span class="line">            count += <span class="number">1</span></span><br><span class="line">            primes[index] = <span class="number">0</span></span><br><span class="line">            <span class="keyword">if</span> i%res[j] == <span class="number">0</span>:</span><br><span class="line">                <span class="keyword">break</span></span><br><span class="line">            <span class="keyword">else</span>:</span><br><span class="line">                j += <span class="number">1</span></span><br><span class="line">                index = <span class="built_in">int</span>(res[j] * i)</span><br><span class="line">    <span class="keyword">return</span> count</span><br></pre></td></tr></table></figure>
<p>观察发现，Euler 方法中，每次循环多出大数取模运算和乘法运算。于是简单测量了一下时长：</p>
<figure class="highlight python"><table><tr><td class="gutter"><pre><span class="line">1</span><br><span class="line">2</span><br><span class="line">3</span><br><span class="line">4</span><br><span class="line">5</span><br><span class="line">6</span><br><span class="line">7</span><br><span class="line">8</span><br><span class="line">9</span><br><span class="line">10</span><br><span class="line">11</span><br><span class="line">12</span><br><span class="line">13</span><br><span class="line">14</span><br><span class="line">15</span><br><span class="line">16</span><br><span class="line">17</span><br></pre></td><td class="code"><pre><span class="line">s = []</span><br><span class="line"><span class="keyword">for</span> _ <span class="keyword">in</span> <span class="built_in">range</span>(<span class="number">10</span>):</span><br><span class="line">    lst = [<span class="number">0</span> <span class="keyword">for</span> _ <span class="keyword">in</span> <span class="built_in">range</span>(<span class="number">10000000</span>)]</span><br><span class="line">    t0 = time.perf_counter()</span><br><span class="line">    lst[<span class="number">5000000</span>] = <span class="number">1</span></span><br><span class="line">    t1 = time.perf_counter()</span><br><span class="line">    dt1 = t1 - t0</span><br><span class="line">    lst = [<span class="number">0</span> <span class="keyword">for</span> _ <span class="keyword">in</span> <span class="built_in">range</span>(<span class="number">10000000</span>)]</span><br><span class="line">    t0 = time.perf_counter()</span><br><span class="line">    lst[<span class="number">5000000</span>] = <span class="number">1</span></span><br><span class="line">    a = <span class="number">399999</span> % <span class="number">19</span></span><br><span class="line">    b = <span class="number">19</span> * <span class="number">399999</span></span><br><span class="line">    t1 = time.perf_counter()</span><br><span class="line">    dt2 = t1 - t0</span><br><span class="line">    ic(<span class="string">f&#x27;<span class="subst">&#123;dt2/dt1:<span class="number">.2</span>f&#125;</span>&#x27;</span>)</span><br><span class="line">    s.append(dt2/dt1)</span><br><span class="line">ic(<span class="built_in">sum</span>(s)/<span class="built_in">len</span>(s))</span><br></pre></td></tr></table></figure>
<p>执行多次以后发现，时长比有着一定范围的波动，最高可到 2 倍以上：</p>
<figure class="highlight bash"><table><tr><td class="gutter"><pre><span class="line">1</span><br><span class="line">2</span><br><span class="line">3</span><br><span class="line">4</span><br><span class="line">5</span><br><span class="line">6</span><br><span class="line">7</span><br><span class="line">8</span><br><span class="line">9</span><br><span class="line">10</span><br><span class="line">11</span><br><span class="line">12</span><br><span class="line">13</span><br><span class="line">14</span><br><span class="line">15</span><br><span class="line">16</span><br><span class="line">17</span><br><span class="line">18</span><br><span class="line">19</span><br><span class="line">20</span><br><span class="line">21</span><br><span class="line">22</span><br><span class="line">23</span><br><span class="line">24</span><br></pre></td><td class="code"><pre><span class="line">ic| f<span class="string">&#x27;&#123;dt2/dt1:.2f&#125;&#x27;</span>: <span class="string">&#x27;2.20&#x27;</span></span><br><span class="line">ic| f<span class="string">&#x27;&#123;dt2/dt1:.2f&#125;&#x27;</span>: <span class="string">&#x27;1.37&#x27;</span></span><br><span class="line">ic| f<span class="string">&#x27;&#123;dt2/dt1:.2f&#125;&#x27;</span>: <span class="string">&#x27;1.44&#x27;</span></span><br><span class="line">ic| f<span class="string">&#x27;&#123;dt2/dt1:.2f&#125;&#x27;</span>: <span class="string">&#x27;1.09&#x27;</span></span><br><span class="line">ic| f<span class="string">&#x27;&#123;dt2/dt1:.2f&#125;&#x27;</span>: <span class="string">&#x27;1.32&#x27;</span></span><br><span class="line">ic| f<span class="string">&#x27;&#123;dt2/dt1:.2f&#125;&#x27;</span>: <span class="string">&#x27;1.16&#x27;</span></span><br><span class="line">ic| f<span class="string">&#x27;&#123;dt2/dt1:.2f&#125;&#x27;</span>: <span class="string">&#x27;1.33&#x27;</span></span><br><span class="line">ic| f<span class="string">&#x27;&#123;dt2/dt1:.2f&#125;&#x27;</span>: <span class="string">&#x27;1.22&#x27;</span></span><br><span class="line">ic| f<span class="string">&#x27;&#123;dt2/dt1:.2f&#125;&#x27;</span>: <span class="string">&#x27;1.33&#x27;</span></span><br><span class="line">ic| f<span class="string">&#x27;&#123;dt2/dt1:.2f&#125;&#x27;</span>: <span class="string">&#x27;1.18&#x27;</span></span><br><span class="line">ic| sum(s)/len(s): 1.3642452780259806</span><br><span class="line"></span><br><span class="line">ic| f<span class="string">&#x27;&#123;dt2/dt1:.2f&#125;&#x27;</span>: <span class="string">&#x27;1.46&#x27;</span></span><br><span class="line">ic| f<span class="string">&#x27;&#123;dt2/dt1:.2f&#125;&#x27;</span>: <span class="string">&#x27;1.20&#x27;</span></span><br><span class="line">ic| f<span class="string">&#x27;&#123;dt2/dt1:.2f&#125;&#x27;</span>: <span class="string">&#x27;1.18&#x27;</span></span><br><span class="line">ic| f<span class="string">&#x27;&#123;dt2/dt1:.2f&#125;&#x27;</span>: <span class="string">&#x27;1.50&#x27;</span></span><br><span class="line">ic| f<span class="string">&#x27;&#123;dt2/dt1:.2f&#125;&#x27;</span>: <span class="string">&#x27;1.00&#x27;</span></span><br><span class="line">ic| f<span class="string">&#x27;&#123;dt2/dt1:.2f&#125;&#x27;</span>: <span class="string">&#x27;1.69&#x27;</span></span><br><span class="line">ic| f<span class="string">&#x27;&#123;dt2/dt1:.2f&#125;&#x27;</span>: <span class="string">&#x27;1.33&#x27;</span></span><br><span class="line">ic| f<span class="string">&#x27;&#123;dt2/dt1:.2f&#125;&#x27;</span>: <span class="string">&#x27;1.33&#x27;</span></span><br><span class="line">ic| f<span class="string">&#x27;&#123;dt2/dt1:.2f&#125;&#x27;</span>: <span class="string">&#x27;1.62&#x27;</span></span><br><span class="line">ic| f<span class="string">&#x27;&#123;dt2/dt1:.2f&#125;&#x27;</span>: <span class="string">&#x27;1.18&#x27;</span></span><br><span class="line">ic| sum(s)/len(s): 1.3493646302728954</span><br><span class="line"></span><br></pre></td></tr></table></figure>
<p>平均水平在 1.3 左右。</p>
<p>所以采用了 Numpy 中的取模运算：</p>
<figure class="highlight python"><table><tr><td class="gutter"><pre><span class="line">1</span><br><span class="line">2</span><br><span class="line">3</span><br><span class="line">4</span><br><span class="line">5</span><br><span class="line">6</span><br><span class="line">7</span><br><span class="line">8</span><br><span class="line">9</span><br><span class="line">10</span><br><span class="line">11</span><br><span class="line">12</span><br><span class="line">13</span><br><span class="line">14</span><br><span class="line">15</span><br><span class="line">16</span><br><span class="line">17</span><br><span class="line">18</span><br><span class="line">19</span><br><span class="line">20</span><br><span class="line">21</span><br><span class="line">22</span><br></pre></td><td class="code"><pre><span class="line"><span class="function"><span class="keyword">def</span> <span class="title">Euler_improved_numpy_with_mod</span>(<span class="params">n</span>):</span></span><br><span class="line">    primes = np.ones(n + <span class="number">1</span>)</span><br><span class="line">    primes[<span class="number">0</span>] = primes[<span class="number">1</span>] = <span class="number">0</span></span><br><span class="line">    res = np.zeros(<span class="built_in">int</span>(n / <span class="number">2</span>))</span><br><span class="line">    count = <span class="number">0</span></span><br><span class="line">    c = <span class="number">0</span></span><br><span class="line">    <span class="keyword">for</span> i <span class="keyword">in</span> <span class="built_in">range</span>(<span class="number">2</span>, n + <span class="number">1</span>):</span><br><span class="line">        <span class="keyword">if</span> primes[i]:</span><br><span class="line">            res[c] = i</span><br><span class="line">            c += <span class="number">1</span></span><br><span class="line">        </span><br><span class="line">        j = <span class="number">0</span></span><br><span class="line">        index = <span class="built_in">int</span>(res[j] * i)</span><br><span class="line">        <span class="keyword">while</span> j &lt; c <span class="keyword">and</span> index &lt;= n:</span><br><span class="line">            count += <span class="number">1</span></span><br><span class="line">            primes[index] = <span class="number">0</span></span><br><span class="line">            <span class="keyword">if</span> np.mod(i,res[j]) == <span class="number">0</span>:</span><br><span class="line">                <span class="keyword">break</span></span><br><span class="line">            <span class="keyword">else</span>:</span><br><span class="line">                j += <span class="number">1</span></span><br><span class="line">                index = <span class="built_in">int</span>(res[j] * i)</span><br><span class="line">    <span class="keyword">return</span> count</span><br></pre></td></tr></table></figure>
<p>我从 StackOverflow 上搜索有关问题，又对埃氏筛进行了改进，同时发现有大佬给出了几种超快的算法<sup id="fnref:3"><a href="#fn:3" rel="footnote">3</a></sup>，我选取了其中两种，一种是纯 Python 的，另一种是使用 Numpy 的：</p>
<figure class="highlight python"><table><tr><td class="gutter"><pre><span class="line">1</span><br><span class="line">2</span><br><span class="line">3</span><br><span class="line">4</span><br><span class="line">5</span><br><span class="line">6</span><br><span class="line">7</span><br><span class="line">8</span><br><span class="line">9</span><br><span class="line">10</span><br><span class="line">11</span><br></pre></td><td class="code"><pre><span class="line"><span class="function"><span class="keyword">def</span> <span class="title">Eratosthenes_numpy_improved</span>(<span class="params">n</span>):</span></span><br><span class="line">    primes = np.ones(n + <span class="number">1</span>)</span><br><span class="line">    primes[<span class="number">0</span>] = primes[<span class="number">1</span>] = <span class="number">0</span></span><br><span class="line">    <span class="keyword">for</span> i <span class="keyword">in</span> <span class="built_in">range</span>(<span class="number">2</span>, <span class="built_in">int</span>(n ** <span class="number">0.5</span>) + <span class="number">1</span>):</span><br><span class="line">        <span class="keyword">if</span> primes[i]:</span><br><span class="line">            primes[i * i::i] = <span class="number">0</span></span><br><span class="line">    c = <span class="number">0</span></span><br><span class="line">    <span class="keyword">for</span> i <span class="keyword">in</span> primes:</span><br><span class="line">        <span class="keyword">if</span> i:</span><br><span class="line">            c += <span class="number">1</span></span><br><span class="line">    <span class="keyword">return</span> c</span><br></pre></td></tr></table></figure>
<p>使用 Numpy 的切片运算可以显著提高速度，Stackoverflow 上的方法也用到了切片运算：</p>
<figure class="highlight python"><table><tr><td class="gutter"><pre><span class="line">1</span><br><span class="line">2</span><br><span class="line">3</span><br><span class="line">4</span><br><span class="line">5</span><br><span class="line">6</span><br><span class="line">7</span><br><span class="line">8</span><br><span class="line">9</span><br><span class="line">10</span><br><span class="line">11</span><br><span class="line">12</span><br><span class="line">13</span><br><span class="line">14</span><br><span class="line">15</span><br><span class="line">16</span><br><span class="line">17</span><br><span class="line">18</span><br><span class="line">19</span><br></pre></td><td class="code"><pre><span class="line"><span class="function"><span class="keyword">def</span> <span class="title">primes</span>(<span class="params">n</span>):</span></span><br><span class="line">    sieve = [<span class="literal">True</span>] * n</span><br><span class="line">    <span class="keyword">for</span> i <span class="keyword">in</span> <span class="built_in">range</span>(<span class="number">3</span>,<span class="built_in">int</span>(n**<span class="number">0.5</span>)+<span class="number">1</span>,<span class="number">2</span>):</span><br><span class="line">        <span class="keyword">if</span> sieve[i]:</span><br><span class="line">            sieve[i*i::<span class="number">2</span>*i]=[<span class="literal">False</span>]*((n-i*i-<span class="number">1</span>)//(<span class="number">2</span>*i)+<span class="number">1</span>)</span><br><span class="line">    res =  [<span class="number">2</span>] + [i <span class="keyword">for</span> i <span class="keyword">in</span> <span class="built_in">range</span>(<span class="number">3</span>,n,<span class="number">2</span>) <span class="keyword">if</span> sieve[i]]</span><br><span class="line">    <span class="keyword">return</span> <span class="built_in">len</span>(res)</span><br><span class="line">    </span><br><span class="line"><span class="function"><span class="keyword">def</span> <span class="title">primesfrom2to</span>(<span class="params">n</span>):</span></span><br><span class="line">    sieve = np.ones(n//<span class="number">3</span> + (n%<span class="number">6</span>==<span class="number">2</span>), </span><br><span class="line">                    dtype=<span class="built_in">bool</span>)</span><br><span class="line">    <span class="keyword">for</span> i <span class="keyword">in</span> <span class="built_in">range</span>(<span class="number">1</span>,<span class="built_in">int</span>(n**<span class="number">0.5</span>)//<span class="number">3</span>+<span class="number">1</span>):</span><br><span class="line">        <span class="keyword">if</span> sieve[i]:</span><br><span class="line">            k=<span class="number">3</span>*i+<span class="number">1</span>|<span class="number">1</span></span><br><span class="line">            sieve[k*k//<span class="number">3</span>::<span class="number">2</span>*k] = <span class="literal">False</span></span><br><span class="line">            sieve[k*(k-<span class="number">2</span>*(i&amp;<span class="number">1</span>)+<span class="number">4</span>)//<span class="number">3</span>::<span class="number">2</span>*k] = <span class="literal">False</span></span><br><span class="line">    res =  np.r_[<span class="number">2</span>,<span class="number">3</span>,((<span class="number">3</span>*np.nonzero(sieve)[<span class="number">0</span>][<span class="number">1</span>:]+<span class="number">1</span>)|<span class="number">1</span>)]</span><br><span class="line">    <span class="keyword">return</span> <span class="built_in">len</span>(res)</span><br><span class="line"></span><br></pre></td></tr></table></figure>
<p>使用不同的数量级对这些函数的速度进行测试，发现大多数函数运行时间基本是 10 倍变化，这是合理的:</p>
<figure class="highlight python"><table><tr><td class="gutter"><pre><span class="line">1</span><br><span class="line">2</span><br><span class="line">3</span><br><span class="line">4</span><br><span class="line">5</span><br><span class="line">6</span><br><span class="line">7</span><br></pre></td><td class="code"><pre><span class="line"><span class="comment"># 测试函数</span></span><br><span class="line"><span class="function"><span class="keyword">def</span> <span class="title">prime_list</span>(<span class="params">func, n</span>):</span></span><br><span class="line">    t0 = time.perf_counter()</span><br><span class="line">    res = func(n)</span><br><span class="line">    t1 = time.perf_counter()</span><br><span class="line">    dt = t1 - t0</span><br><span class="line">    <span class="keyword">return</span> dt, res</span><br></pre></td></tr></table></figure>
<p>时间结果附在文后<sup id="fnref:4"><a href="#fn:4" rel="footnote">4</a></sup>。</p>
<figure class="highlight python"><table><tr><td class="gutter"><pre><span class="line">1</span><br><span class="line">2</span><br><span class="line">3</span><br><span class="line">4</span><br><span class="line">5</span><br><span class="line">6</span><br><span class="line">7</span><br><span class="line">8</span><br><span class="line">9</span><br><span class="line">10</span><br><span class="line">11</span><br><span class="line">12</span><br><span class="line">13</span><br><span class="line">14</span><br><span class="line">15</span><br><span class="line">16</span><br><span class="line">17</span><br><span class="line">18</span><br><span class="line">19</span><br><span class="line">20</span><br><span class="line">21</span><br><span class="line">22</span><br><span class="line">23</span><br><span class="line">24</span><br><span class="line">25</span><br><span class="line">26</span><br><span class="line">27</span><br><span class="line">28</span><br></pre></td><td class="code"><pre><span class="line"><span class="comment"># 绘图代码</span></span><br><span class="line"><span class="keyword">import</span> pandas <span class="keyword">as</span> pd</span><br><span class="line"><span class="keyword">import</span> matplotlib.pyplot <span class="keyword">as</span> plt</span><br><span class="line"><span class="keyword">from</span> icecream <span class="keyword">import</span> ic</span><br><span class="line"><span class="keyword">import</span> numpy <span class="keyword">as</span> np</span><br><span class="line"></span><br><span class="line"><span class="function"><span class="keyword">def</span> <span class="title">to_float</span>(<span class="params">series</span>):</span></span><br><span class="line">    series_new = pd.Series(dtype=<span class="string">&#x27;float64&#x27;</span>)</span><br><span class="line">    <span class="keyword">for</span> i <span class="keyword">in</span> <span class="built_in">range</span>(<span class="built_in">len</span>(series)):</span><br><span class="line">        series_new[<span class="built_in">str</span>(i)] = np.log(<span class="built_in">float</span>(series.iloc[i]))+<span class="number">14</span></span><br><span class="line">    argmax = <span class="built_in">int</span>(series_new.idxmax())</span><br><span class="line">    argmin = <span class="built_in">int</span>(series_new.idxmin())</span><br><span class="line">    <span class="keyword">return</span> &#123;<span class="string">&#x27;y&#x27;</span>:series_new,<span class="string">&#x27;idxmax&#x27;</span>:argmax,<span class="string">&#x27;max&#x27;</span>:series_new[argmax],</span><br><span class="line">            <span class="string">&#x27;idxmin&#x27;</span>:argmin,<span class="string">&#x27;min&#x27;</span>:series_new[argmin]&#125;</span><br><span class="line"></span><br><span class="line">df = pd.read_csv(<span class="string">&#x27; 素数筛速度统计.csv&#x27;</span>,dtype = <span class="string">&#x27;object&#x27;</span>).\\</span><br><span class="line">                  drop(columns=[<span class="string">&#x27;numba(jit no python)&#x27;</span>],axis=<span class="number">1</span>)</span><br><span class="line">nojit = df[<span class="number">0</span>:<span class="number">11</span>]</span><br><span class="line">funcs = df[<span class="string">&#x27;function&#x27;</span>][<span class="number">0</span>:<span class="number">11</span>]</span><br><span class="line">plt.figure(figsize=(<span class="number">10</span>,<span class="number">10</span>))</span><br><span class="line"><span class="keyword">for</span> i <span class="keyword">in</span> <span class="built_in">range</span>(<span class="number">8</span>,<span class="number">0</span>,-<span class="number">1</span>):</span><br><span class="line">    res = to_float(nojit.loc[:,<span class="built_in">str</span>(<span class="number">10</span> ** i)])</span><br><span class="line">    plt.barh(funcs,res[<span class="string">&#x27;y&#x27;</span>])</span><br><span class="line">plt.yticks(fontsize = <span class="number">10</span>)</span><br><span class="line">plt.subplots_adjust(left=<span class="number">0.25</span>)</span><br><span class="line">plt.legend(loc = <span class="string">&#x27;upper right&#x27;</span>,labels = [<span class="number">10</span> ** i <span class="keyword">for</span> i <span class="keyword">in</span> <span class="built_in">range</span>(<span class="number">8</span>,<span class="number">0</span>,-<span class="number">1</span>)])</span><br><span class="line">plt.savefig(<span class="string">&#x27;nojit.svg&#x27;</span>)</span><br><span class="line">plt.show()</span><br></pre></td></tr></table></figure>
<p>得到图像：</p>
<div class="tag-plugin image"><div class="image-bg"><img class="lazy" src="" data-src="nojit.svg"/><a class="image-download blur" style="opacity:0" href="nojit.svg"><svg class="icon" style="width: 1em; height: 1em;vertical-align: middle;fill: currentColor;overflow: hidden;" viewBox="0 0 1024 1024" version="1.1" xmlns="http://www.w3.org/2000/svg" p-id="3734"><path d="M561.00682908 685.55838913a111.03077546 111.03077546 0 0 1-106.8895062 0L256.23182837 487.72885783a55.96309219 55.96309219 0 0 1 79.13181253-79.18777574L450.70357448 523.88101491V181.55477937a55.96309219 55.96309219 0 0 1 111.92618438 0v344.06109173l117.07478902-117.07478901a55.96309219 55.96309219 0 0 1 79.13181252 79.18777574zM282.81429711 797.1487951h447.70473912a55.96309219 55.96309219 0 0 1 0 111.92618438H282.81429711a55.96309219 55.96309219 0 0 1 0-111.92618438z" p-id="3735"></path></svg></a></div></div>
<p>这里对运行时间进行了处理，来更方便地展示：</p>
<p class='katex-block'><span class="katex-display"><span class="katex"><span class="katex-mathml"><math><semantics><mrow><msup><mi>t</mi><mo mathvariant="normal">′</mo></msup><mo>=</mo><mi>ln</mi><mo>⁡</mo><mi>t</mi><mo>+</mo><mn>14</mn></mrow><annotation encoding="application/x-tex">t&#x27;=\ln t+14
</annotation></semantics></math></span><span class="katex-html" aria-hidden="true"><span class="base"><span class="strut" style="height:0.801892em;vertical-align:0em;"></span><span class="mord"><span class="mord mathdefault">t</span><span class="msupsub"><span class="vlist-t"><span class="vlist-r"><span class="vlist" style="height:0.801892em;"><span style="top:-3.113em;margin-right:0.05em;"><span class="pstrut" style="height:2.7em;"></span><span class="sizing reset-size6 size3 mtight"><span class="mord mtight"><span class="mord mtight">′</span></span></span></span></span></span></span></span></span><span class="mspace" style="margin-right:0.2777777777777778em;"></span><span class="mrel">=</span><span class="mspace" style="margin-right:0.2777777777777778em;"></span></span><span class="base"><span class="strut" style="height:0.77777em;vertical-align:-0.08333em;"></span><span class="mop">ln</span><span class="mspace" style="margin-right:0.16666666666666666em;"></span><span class="mord mathdefault">t</span><span class="mspace" style="margin-right:0.2222222222222222em;"></span><span class="mbin">+</span><span class="mspace" style="margin-right:0.2222222222222222em;"></span></span><span class="base"><span class="strut" style="height:0.64444em;vertical-align:0em;"></span><span class="mord">1</span><span class="mord">4</span></span></span></span></span></p>
<p>可以看到，在参考 StackOverflow 前的方法中，仅仅是使用 Numpy 生成数组并不能加快速度，反而会变得很慢。同时，虽然 Euler 方法减少了重复，但实际运行速度并没有加快，可能有循环步骤中操作过多的原因。而且在小规模时，使用 Numpy 的 Eratosthenes numpy improved 和 primesfrom2to 并没有表现出优势，大规模时优势才凸显出来，而使用原生列表的 primes 始终保持着较好的速度，这可能是 Numpy 的原因。</p>
<p>在使用 Numba 进行加速以后，得到新的运行时间：</p>
<figure class="highlight python"><table><tr><td class="gutter"><pre><span class="line">1</span><br><span class="line">2</span><br></pre></td><td class="code"><pre><span class="line"><span class="comment"># 在每个函数前添加</span></span><br><span class="line"><span class="meta"> @jit(<span class="params">nopython = <span class="literal">True</span></span>)</span></span><br></pre></td></tr></table></figure>
<p>需要说明的是，由于 primesfrom2to 无法完全脱离 Python 解释器运行，所以使用了：</p>
<figure class="highlight python"><table><tr><td class="gutter"><pre><span class="line">1</span><br></pre></td><td class="code"><pre><span class="line"><span class="meta">@jit(<span class="params">）</span></span></span><br></pre></td></tr></table></figure>
<p>使用同样的方法进行测试并绘制图像：</p>
<figure class="highlight python"><table><tr><td class="gutter"><pre><span class="line">1</span><br><span class="line">2</span><br><span class="line">3</span><br><span class="line">4</span><br><span class="line">5</span><br><span class="line">6</span><br><span class="line">7</span><br><span class="line">8</span><br><span class="line">9</span><br><span class="line">10</span><br><span class="line">11</span><br><span class="line">12</span><br><span class="line">13</span><br></pre></td><td class="code"><pre><span class="line">jit = df[<span class="number">11</span>:]</span><br><span class="line">plt.figure(figsize=(<span class="number">10</span>,<span class="number">10</span>))</span><br><span class="line"><span class="keyword">for</span> i <span class="keyword">in</span> <span class="built_in">range</span>(<span class="number">8</span>,<span class="number">0</span>,-<span class="number">1</span>):</span><br><span class="line">    res = to_float(jit.loc[:,<span class="built_in">str</span>(<span class="number">10</span> ** i)])</span><br><span class="line">    <span class="keyword">if</span> i == <span class="number">1</span>:</span><br><span class="line">        plt.barh(funcs,res[<span class="string">&#x27;y&#x27;</span>],alpha = <span class="number">0.3</span>)</span><br><span class="line">    <span class="keyword">else</span>:</span><br><span class="line">        plt.barh(funcs,res[<span class="string">&#x27;y&#x27;</span>])</span><br><span class="line">plt.yticks(fontsize = <span class="number">10</span>)</span><br><span class="line">plt.subplots_adjust(left=<span class="number">0.25</span>)</span><br><span class="line">plt.legend(loc = <span class="string">&#x27;upper right&#x27;</span>,labels = [<span class="number">10</span> ** i <span class="keyword">for</span> i <span class="keyword">in</span> <span class="built_in">range</span>(<span class="number">8</span>,<span class="number">0</span>,-<span class="number">1</span>)])</span><br><span class="line">plt.savefig(<span class="string">&#x27;jit.svg&#x27;</span>)</span><br><span class="line">plt.show() </span><br></pre></td></tr></table></figure>
<p>这是由于在规模为 10 时，速度较慢，呈现在图中会严重覆盖其他内容，所以调整了它的透明度，图像中变暗的地方时规模为 10 的情况：</p>
<div class="tag-plugin image"><div class="image-bg"><img class="lazy" src="" data-src="jit.svg"/><a class="image-download blur" style="opacity:0" href="jit.svg"><svg class="icon" style="width: 1em; height: 1em;vertical-align: middle;fill: currentColor;overflow: hidden;" viewBox="0 0 1024 1024" version="1.1" xmlns="http://www.w3.org/2000/svg" p-id="3734"><path d="M561.00682908 685.55838913a111.03077546 111.03077546 0 0 1-106.8895062 0L256.23182837 487.72885783a55.96309219 55.96309219 0 0 1 79.13181253-79.18777574L450.70357448 523.88101491V181.55477937a55.96309219 55.96309219 0 0 1 111.92618438 0v344.06109173l117.07478902-117.07478901a55.96309219 55.96309219 0 0 1 79.13181252 79.18777574zM282.81429711 797.1487951h447.70473912a55.96309219 55.96309219 0 0 1 0 111.92618438H282.81429711a55.96309219 55.96309219 0 0 1 0-111.92618438z" p-id="3735"></path></svg></a></div></div>
<p>可以看到，Numba 的加速效果是显而易见的，速度最慢是 16 左右。由于 primesfrom2to 基本没有变化，而且本身速度很快，所以这里不考虑。考虑其他函数，发现加速后 Numpy 的优势表现出来，而且经过加速，Euler 方法在大规模问题上优于 Eratosthenes 方法。有趣的是，加速以后 Euler 的 improved 方法表现变差，可能是 Numba 对于 <code>append</code> 操作有着较好的优化的原因，比我自己简单的改进更高效。</p>
<p>总的来看，使用 Eratosthenes numpy improved 或者 Numba 加速后的 Euler numpy 就可以在 1 亿规模获得简单好用的效果。</p>
<div id="footnotes"><hr><div id="footnotelist"><ol style="list-style:none; padding-left: 0;"><li id="fn:1"><span style="display: inline-block; vertical-align: top; padding-right: 10px;">1.</span><span style="display: inline-block; vertical-align: top;"><a target="_blank" rel="noopener" href="https://blog.csdn.net/mxYlulu/article/details/81806980">快速筛素数（埃式筛 + 线性筛 +Miller_Rabin 算法）</a></span><a href="#fnref:1" rev="footnote"> ↩</a></li><li id="fn:2"><span style="display: inline-block; vertical-align: top; padding-right: 10px;">2.</span><span style="display: inline-block; vertical-align: top;"><a target="_blank" rel="noopener" href="https://blog.csdn.net/qq_41117236/article/details/81152055">判断素数（一般筛到线性筛）</a></span><a href="#fnref:2" rev="footnote"> ↩</a></li><li id="fn:3"><span style="display: inline-block; vertical-align: top; padding-right: 10px;">3.</span><span style="display: inline-block; vertical-align: top;"><a target="_blank" rel="noopener" href="https://stackoverflow.com/questions/2068372/fastest-way-to-list-all-primes-below-n">Fastest way to list all primes below N</a></span><a href="#fnref:3" rev="footnote"> ↩</a></li><li id="fn:4"><span style="display: inline-block; vertical-align: top; padding-right: 10px;">4.</span><span style="display: inline-block; vertical-align: top;"><a target="_blank" rel="noopener" href="https://gitee.com/stassenger/stassenger/blob/file/%E7%B4%A0%E6%95%B0%E7%AD%9B%E9%80%9F%E5%BA%A6%E7%BB%9F%E8%AE%A1.csv">素数筛速度统计.csv</a></span><a href="#fnref:4" rev="footnote"> ↩</a></li></ol></div></div>

<div class="article-footer fs14"><section id="license"><div class="header"><span>许可协议</span></div><div class="body"><p>本文采用 <a target="_blank" rel="noopener" href="https://creativecommons.org/licenses/by-nc-sa/4.0/">署名-非商业性使用-相同方式共享 4.0 国际</a> 许可协议，转载请注明出处。</p>
</div></section><section id="share"><div class="header"><span>分享文章</span></div><div class="body"><div class="link"><input class="copy-area" readonly="true" id="copy-link" value="https://stassenger.gitee.io/2021/07/08/a%20interesting%20problem%20when%20I%E2%80%98m%20learning%20Euler%20and%20Eratosthenes%20method/" /></div><div class="social-wrap dis-select"><a class="social share-item wechat" onclick="util.toggle(&quot;qrcode-wechat&quot)"><img class="lazy" src="" data-src="https://cdn.jsdelivr.net/gh/cdn-x/placeholder@1.0.1/social/b32ef3da1162a.svg"/></a><a class="social share-item email" href="mailto:?subject=使用 Python 实现素数筛时发现的有趣问题 - Debug Myself&amp;body=https://stassenger.gitee.io/2021/07/08/a%20interesting%20problem%20when%20I%E2%80%98m%20learning%20Euler%20and%20Eratosthenes%20method/"><img class="lazy" src="" data-src="https://cdn.jsdelivr.net/gh/cdn-x/placeholder@1.0.1/social/a1b00e20f425d.svg"/></a><a class="social share-item link" onclick="util.copy(&quot;copy-link&quot;, &quot;复制成功&quot;)"><img class="lazy" src="" data-src="https://cdn.jsdelivr.net/gh/cdn-x/placeholder@1.0.1/social/8411ed322ced6.svg"/></a></div><div class="qrcode" id="qrcode-wechat" style="visibility:hidden;height:0"><img class="lazy" src="" data-src="https://api.qrserver.com/v1/create-qr-code/?size=256x256&data=https://stassenger.gitee.io/2021/07/08/a%20interesting%20problem%20when%20I%E2%80%98m%20learning%20Euler%20and%20Eratosthenes%20method/index.html"/></div></div></section></div>

</article>

<div class="related-wrap" id="read-next"><section class="header cap theme"><span>接下来阅读</span></section><section class="body"><div class="post-title h2"><a href="/2021/06/14/some%20of%20my%20experiences%20about%20cloud%20services/">我使用云服务的一些经历</a></div><div class="post-title fs14"><a href="/2021/07/12/sort-algorithm-reminder-3/">上一篇：排序算法备忘录（三）</a></div></section></div>






  <div class='related-wrap md' id="comments">
    <div class='cmt-title cap theme'>
      快来参与讨论吧
    </div>
    <div class='cmt-body valine'>
      

<div id="valine_container" class="valine_thread"><svg class="loading" style="vertical-align: middle;fill: currentColor;overflow: hidden;" viewBox="0 0 1024 1024" version="1.1" xmlns="http://www.w3.org/2000/svg" p-id="2709"><path d="M832 512c0-176-144-320-320-320V128c211.2 0 384 172.8 384 384h-64zM192 512c0 176 144 320 320 320v64C300.8 896 128 723.2 128 512h64z" p-id="2710"></path></svg></div>

    </div>
  </div>



      
<footer class="page-footer fs12"><hr><div><p>本博客所有文章除特别声明外，均采用 <a target="_blank" rel="noopener" href="https://creativecommons.org/licenses/by-nc-sa/4.0/">CC BY-NC-SA 4.0</a> 许可协议，转载请注明出处。</p>
<p>本站由 <a href="https://stassenger.gitee.io/">@Asante Stassenger</a> 创建，使用 <a target="_blank" rel="noopener" href="https://github.com/xaoxuu/hexo-theme-stellar/tree/1.1.0">Stellar</a> 作为主题。</p>
</div></footer>


      <div class='float-panel mobile-only' style='display:none'>
  <button type='button' class='sidebar-toggle mobile' onclick='sidebar.toggle()'>
    <svg class="icon" viewBox="0 0 1228 1024" version="1.1" xmlns="http://www.w3.org/2000/svg" p-id="2849"><path d="M0 0m97.390199 0l973.901995 0q97.390199 0 97.390199 97.390199l0 0q0 97.390199-97.390199 97.3902l-973.901995 0q-97.390199 0-97.390199-97.3902l0 0q0-97.390199 97.390199-97.390199Z" p-id="2850"></path><path d="M0 389.560798m97.390199 0l973.901995 0q97.390199 0 97.390199 97.390199l0 0q0 97.390199-97.390199 97.3902l-973.901995 0q-97.390199 0-97.390199-97.3902l0 0q0-97.390199 97.390199-97.390199Z" p-id="2851"></path><path d="M0 779.121596m97.390199 0l973.901995 0q97.390199 0 97.390199 97.390199l0 0q0 97.390199-97.390199 97.390199l-973.901995 0q-97.390199 0-97.390199-97.390199l0 0q0-97.390199 97.390199-97.390199Z" p-id="2852"></path></svg>
  </button>
</div>

    </div>
  </div>
  <div class='scripts'>
    <script type="text/javascript">
  stellar = {
    // 懒加载 css https://github.com/filamentgroup/loadCSS
    loadCSS: (href, before, media, attributes) => {
      var doc = window.document;
      var ss = doc.createElement("link");
      var ref;
      if (before) {
        ref = before;
      } else {
        var refs = (doc.body || doc.getElementsByTagName("head")[0]).childNodes;
        ref = refs[refs.length - 1];
      }
      var sheets = doc.styleSheets;
      if (attributes) {
        for (var attributeName in attributes) {
          if (attributes.hasOwnProperty(attributeName)) {
            ss.setAttribute(attributeName, attributes[attributeName]);
          }
        }
      }
      ss.rel = "stylesheet";
      ss.href = href;
      ss.media = "only x";
      function ready(cb) {
        if (doc.body) {
          return cb();
        }
        setTimeout(function () {
          ready(cb);
        });
      }
      ready(function () {
        ref.parentNode.insertBefore(ss, before ? ref : ref.nextSibling);
      });
      var onloadcssdefined = function (cb) {
        var resolvedHref = ss.href;
        var i = sheets.length;
        while (i--) {
          if (sheets[i].href === resolvedHref) {
            return cb();
          }
        }
        setTimeout(function () {
          onloadcssdefined(cb);
        });
      };
      function loadCB() {
        if (ss.addEventListener) {
          ss.removeEventListener("load", loadCB);
        }
        ss.media = media || "all";
      }
      if (ss.addEventListener) {
        ss.addEventListener("load", loadCB);
      }
      ss.onloadcssdefined = onloadcssdefined;
      onloadcssdefined(loadCB);
      return ss;
    },

    // 从 butterfly 和 volantis 获得灵感
    loadScript: (src, opt) => new Promise((resolve, reject) => {
      var script = document.createElement('script')
      script.src = src
      if (opt) {
        for (let key of Object.keys(opt)) {
          script[key] = opt[key]
        }
      } else {
        // 默认异步，如果需要同步，第二个参数传入 {} 即可
        script.async = true
      }
      script.onerror = reject
      script.onload = script.onreadystatechange = function() {
        const loadState = this.readyState
        if (loadState && loadState !== 'loaded' && loadState !== 'complete') return
        script.onload = script.onreadystatechange = null
        resolve()
      }
      document.head.appendChild(script)
    }),

    // https://github.com/jerryc127/hexo-theme-butterfly
    jQuery: (fn) => {
      if (typeof jQuery === 'undefined') {
        stellar.loadScript(stellar.plugins.jQuery).then(fn)
      } else {
        fn()
      }
    }
  };
  stellar.github = 'https://github.com/xaoxuu/hexo-theme-stellar/tree/1.1.0';
  stellar.config = {
    date_suffix: {
      just: '刚刚',
      min: '分钟前',
      hour: '小时前',
      day: '天前',
      month: '个月前',
    },
  };

  // required plugins (only load if needs)
  stellar.plugins = {
    jQuery: 'https://cdn.jsdelivr.net/npm/jquery@3.5.1/dist/jquery.min.js',
    sitesjs: '/js/plugins/sites.js',
    friendsjs: '/js/plugins/friends.js',
  };

  // optional plugins
  if ('true' == 'true') {
    stellar.plugins.lazyload = Object.assign({"enable":true,"js":"https://cdn.jsdelivr.net/npm/vanilla-lazyload@17.3.1/dist/lazyload.min.js","transition":"blur"});
  }
  if ('true' == 'true') {
    stellar.plugins.swiper = Object.assign({"enable":true,"css":"https://unpkg.com/swiper/swiper-bundle.min.css","js":"https://unpkg.com/swiper/swiper-bundle.min.js"});
  }
  if ('' == 'true') {
    stellar.plugins.scrollreveal = Object.assign({"enable":null,"js":"https://cdn.jsdelivr.net/npm/scrollreveal@4.0.9/dist/scrollreveal.min.js","distance":"8px","duration":500,"interval":100,"scale":1});
  }
  if ('true' == 'true') {
    stellar.plugins.preload = Object.assign({"enable":true,"service":"flying_pages","instant_page":"https://cdn.jsdelivr.net/gh/volantis-x/cdn-volantis@4.1.2/js/instant_page.js","flying_pages":"https://cdn.jsdelivr.net/gh/gijo-varghese/flying-pages@2.1.2/flying-pages.min.js"});
  }
</script>

<!-- required -->

  
<script src="/js/main.js" async></script>



<!-- optional -->

  <script>
  function getEmojiMaps() {
    function emoji(path, idx, ext) {
      return path + "/" + path + "-" + idx + "." + ext;
    }
    var emojiMaps = {};
    for (var i = 1; i <= 54; i++) {
      emojiMaps['tieba-' + i] = emoji('tieba', i, 'png');
    }
    for (var i = 1; i <= 101; i++) {
      emojiMaps['qq-' + i] = emoji('qq', i, 'gif');
    }
    for (var i = 1; i <= 116; i++) {
      emojiMaps['aru-' + i] = emoji('aru', i, 'gif');
    }
    for (var i = 1; i <= 125; i++) {
      emojiMaps['twemoji-' + i] = emoji('twemoji', i, 'png');
    }
    for (var i = 1; i <= 4; i++) {
      emojiMaps['weibo-' + i] = emoji('weibo', i, 'png');
    }
    return emojiMaps;
  }
  function load_comment(){
    if(!document.getElementById("valine_container"))return;
    stellar.loadScript('https://cdn.jsdelivr.net/gh/XuxuGood/simple-blog-cdn@main/js/Valine.min.js', {defer:true}).then(function () {
      const el = document.getElementById("valine_container");
      var path = el.getAttribute('comment_id');
      const placeholder = "";
      if (!path) {
        path = decodeURI(window.location.pathname);
      }
      var valine = new Valine();
      valine.init(Object.assign({"js":"https://cdn.jsdelivr.net/gh/XuxuGood/simple-blog-cdn@main/js/Valine.min.js","appId":"Gmd8dpj69RW0YYzTOgchjKDG-9Nh9j0Va","appKey":"8XTx3FrEGvL5fwadDKGc0sHf","placeholder":null,"visitor":true,"requiredFields":["nick","mail"],"enableQQ":true,"recordIP":false,"avatar":"robohash","pageSize":10,"lang":"zh-cn","highlight":true,"mathJax":false,"tagMeta":["博主","小伙伴","访客"],"metaPlaceholder":{"nick":"昵称/QQ号(必填)","mail":"邮箱（必填）"},"master":["6783037F2DF30EAB99F9FC256157D875"],"friends":["6783037F2DF30EAB99F9FC256157D875"]}, {
        el: '#valine_container',
        path: path,
        placeholder: placeholder,
        emojiCDN: 'https://cdn.jsdelivr.net/gh/volantis-x/cdn-emoji/valine/',
        emojiMaps: getEmojiMaps(),
      }));
    });
  }
  window.addEventListener('DOMContentLoaded', (event) => {
    console.log('DOM fully loaded and parsed');
    load_comment();
  });

</script>




<!-- inject -->


  </div>
</body>
</html>
